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ABSTRACT 


This  report  addresses  the  development  of  simulation  tools  which  are  suitable  for  the 
design  and  analysis  of  optimum  imaging  scenarios  for  adaptive  optics  systems.  In  particular,  we 
discuss  the  physical  basis  and  suitability  of  the  modular  CAOS  system  for  such  studies.  A 
desciiption  of  its  essential  features  and  the  basic  philosophy  behind  its  use  is  included  together 
with  an  illustrative  design.  Summary  documentation  is  also  included.  We  conclude  that  CAOS  is 
a  powerful  system  with  considerable  flexibility  and  possessing  considerable  potential  for 
adaptive  optics  study  and  analysis. 

A  second  part  of  the  report  is  given  in  the  form  of  a  paper  submitted  June  8lh  2000  to 
Journal  of  the  Optical  Society  of  America  A  —  “Variational  solution  for  modal  wavefront 
projection  functions  of  minimum  error  norm”.  This  describes  the  derivation  of  vector  functions 
which  are  orthogonal  to  the  gradients  of  the  wavefront  modes.  Their  property  of  minimum  noise 
pi  opagation  allows  an  optimal  method  for  estimating  atmospherically  distorted  wavefronts 
thiough  direct  integration  of  the  wavefront  slope  measurements. 


Both  solutions  have  in  fact  produced  a  variety  of  new  problems.  In 

of  the  artificial  star  (the  so-called  cone  effect  [5])  and  the  bhndnes  <o )  he  rp-  ‘  »f  “ 

wavefront  (lip-tilt  indeterminacy  [6])  are  major  cons, derations.  Th cos, o  such 
considerable  AO  systems  which  will  employ  LGS  must  therefore  be  efficiently  implemented 
Another  topic  of  quite  general  interest  in  AO  relates  to  the  optimal  use  of  the  i  ncoming  photon 
Box  Onefe  rarely  afforded  .he  luxury  of  high  signal-noise  ratios  and  this  implies  an  ongoing 
interest  in  the  study  of  innovative  wavefront  sensing  techniques.  Opt, maOlse  ot The  avada 
nhoton  flux  would  allow  us  to  increase  the  limiting  magnitude  ot  NGS  AO  (inutasin=  tne  .  k 
eoveragel'orto  reduce  the  necessary  laser  power.  For  these  applications  and  others  a  powerful, 
realistic  and  versatile  simulation  tool  is  desiiable. 


CAOS  (Code  for  adaptive  optics  systems)  is  a  purpose-designed  software  tool  for  the  simulation 
LAUS  e  .  Written  in  the  1DL  uro^ramming  language,  there  have  been  a  numbei 

°f  Contributor! 'and  C AOS  iTfreelv  ayailable  to  workers  in  “adaptive  optics.  The  development  ot 
CAOS  came  from  a  need  to  provrde  reliable  and  comprehensive  software  a 

iwo  nrinctaUims  The  firs,  is  to  accnrarcly  simulate  real  AO  sysrems  in  order  ,0  evaluate  and 
characterize  system  performance.  The  second  is  to  allow  the  evaluation  ot  new  techniques. 

designs  and  studies. 


1 .2  Features  of  CAOS 


CAOS  has  two  key  features  which  make  it  easy  to  use  - 

•  A  modular  structure 

•  A  graphical  application  builder 


1.2.1  Modules 

Modules  are  rhe  basic  building  blocks  of  the  simulation.  Each  module  cormsponds  to  a  partien'ar 
identifiable  aspect  of  the  adaptive  optics  system  or  rhe  interpretation  ot  ,1s  results.  Chaos 
version  2.0  has  the  following  modules- 


Atmospheric  building 
Beam  splitter 


ATM 

BSP 


Centroiding  calculus  CEN 

Calibration  fibre  CFB 

Command  sequencer  CSQ 

Data  display  DIS 

Geometrical  propagation  DMI 

interferometric  beam  combiner  IBC 

Imaging  device  IMG 

Laser  definition  LAS 

Make  calibration  data  MCA 

Sodium  layer  NLS 

Wavefront  reconstruction  REC 

Data  saving  SAV 

Shack  Hartmann  sensor  SHS 

Source  definition  SRC 

Structure  function  calculus  STF 

Tip  -tilt  centroiding  TCE 

Time  filtering  TFL 

Tip-tilt  Mirror  TTM 

Wavefront  adding  WFA 


Details  of  some  of  the  most  important  modules  are  given  later  in  the  report. 

1.2.2  Philosophy  of  CAOS 

Each  module  is  defined  by  a  standard  group  of  function  calls,  a  set  of  parameters  which  are 
required  for  the  module's  correct  execution  and  a  set  of  predetermined  inputs  and  outputs. 

The  basic  philosophy  is  that  the  user  not  be  reciuired  to  perform  any  coding  but  need  only 
combine  the  modules  in  a  meaningful  way  and  supply  meaningful  values  of  the  parameters  to 
each  of  the  modules.  In  this  way,  and  within  the  constraints  imposed  by  the  software,  the 
designer  may  concentrate  on  the  concepts.  To  assist  with  this  simulation  design,  a  graphical 
interface  (gui)  termed  the  application  builder  is  supplied. 

1.3  The  application  builder  (AB) 

The  application  builder  is  a  graphical  tool  provided  to  aid  the  design  of  simulation  programs.  It 
used  only  during  the  design  phase.  After  that  the  simulation  program  will  run  as  a  stand-alone 
application  with  no  subsequent  interaction  with  the  AB. 

A  schematic  diagram  of  an  arbitrary  simulation  program  designed  via  the  AB.  together  with  the 
module  list  scroll-down  menu  of  the  AB,  is  shown  in  Fig.  I .  The  simple  interface  consists  of  a 
worksheet  -  essentially  a  collection  of  contiguous  cells  -  with  some  pull-down  menus.  Each  cell 
in  the  worksheet  may  accomodate  a  module. 
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Figure  1:  CAOS  Worksheet  with  pull-down  menus 

Essentially,  modules  are  selected  from  the  pull-down  "Modules"  menu  and  dropped  into  a  chosen 
cell  position  on  the  worksheet.  From  the  point  of  view  of  the  application  builder,  each  module 
can  only  be  executed  if  the  required  set  of  inputs  to  the  module  are  available.  These  inputs  are 
geneially  provided  by  the  outputs  of  other  modules.  Accordingly,  the  input  and  output  ends  of 
each  module  are  colour  coded  to  indicate  compatibility.  The  basic  appearance  of  a  module  as  it 
appears  on  the  worksheet  is  given  below  in  figure  2. 


INPUT 

END 


Figure  2:  Simulation  Module 

The  colour  of  the  input  and  output  simply  specifies  compatibility  between  modules.  These  cannot 
be  changed  or  edited.  Each  output  may  be  joined  to  an  arbitrary  number  of  input  handles  but 


mpirt  handles  may  only  receive  a  single  link  from  an  output  handle.  Conceptually  speaking  the 
ody  of  the  module  contains  the  parameters  which  are  required  for  the  module  to  execute  once 
e  esign  is  complete.  These  can  be  set  and  edited  by  clicking  with  the  mouse  on  the  body  of  the 
module  and  this  invokes  the  corresponding  "parameter  setting"  GUI  -  i.e.  a  graphical  window 
containing  the  current  values  of  the  required  parameters  as  editable  text  fields. 

Links  are  then  defined  between  the  output  of  a  given  module  and  the  input  of  another  thus 
defining  the  data  flow  within  a  simulation  project.  In  the  simple  example  below  in  figure  3  we 
see  the  source  module  (which  defines  the  type  of  luminous  object  we  are  looking  at)  and  the 
atmospheric  builder  module  being  linked  with  the  geometric  propagation  module  In  this 
particular  case,  these  three  modules  would  be  suitable  to  simulate  the  propagation  of  the 

wavefront  from  a  chosen  source  through  an  atmosphere  of  selected  spatiotemporal  behaviour  to 
the  pupil  plane  of  the  telescope. 


Note  that  the  output  colour  of  each  module  in  the  sequence  matches  the  input  colour.  Grey  is 
used  to  indicate  compatibility  with  all  modules.  Cross-hatching  at  the  input  indicates  that  no 
other  modules  may  feed  into  the  given  module. 

The  links  (indicated  in  this  case  simply  by  straight  lines)  are  made  by  clicking  on  the  output  of 
one  module  with  the  mouse  followed  by  clicking  on  the  input  of  the  other. 

After  the  design  is  completed,  the  block  diagram  of  linked  modules  is  analysed  by  the  application 
builder  and  the  IDL  code  implementing  the  simulation  program  is  generated.  The  whole 

simulation  can  be  saved  as  a  project  which  may  be  executed  and  opened  for  later  modification  if 
desired. 


1.4  CAOS  Design  -  basic  summary 

In  summary,  the  simulation  design  procedure  may  be  reduced  to  the  following  sequence  - 
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TASK 

ACTION 

Select  required  modules 

From  pull-down  menu  in  the  application 
builder,  select  and  drop  modules  into  cells  on 
worksheet 

Specify  the  values  of  the  parameters  required 
by  each  module 

Click  on  body  of  selected  module  to  invoke 
GUI  parameter  setting  window.  Set  parameter 
values  and  save 

Form  legal  links  between  modules 

Click  on  the  output  end  of  module  A  followed 
by  a  click  on  the  input  end  of  module  B.  This 
links  A  and  B 

Save  project 

Select  from  pull-down  menu 

1.5  Running  CAOS  simulations 

CAOS  simulations  are  run  as  normal  IDL  projects  from  the  command  line  IDL  interface.  Any 
inconsistencies  in  the  parameter  setting  of  each  module  are  normally  detected  prior  to  this  stage 
by  the  application  builder  which  will  indicate  where  the  inconsistency  has  occurred  and  suggest 
alternative  values. 

1.6  Platform  dependence,  reliability  and  bugs 

CAOS  has  been  primarily  developed  to  run  on  a  UNIX/LINUX  platform.  As  freeware,  the 
developers  make  no  assurances  as  to  reliability  and  existence  of  bugs.  The  reliability  of  CAOS 
linux  v2.0  is  however  considerably  better  than  previous  version.  CAOS  2.0  for  MsWindows95- 
98  exists  but  reliability  is  currently  not  as  good. 

1.7  Obtaining  CAOS 

Scientific  workers  who  wish  to  use  CAOS  should  mail  marcel@arcetri.astrs.it.  On 
receiving  permission,  they  may  then  download  the  software  from  the  corresponding  ftp  site. 


1.8  EXAMPLE  APPLICATION 

In  this  section,  simply  to  illustrate  the  concept  more  clearly,  we  present  a  conceptual  design  and 
result  obtained  from  a  relatively  complex  LGS  adaptive  optics  system.  Full  details  are  available 
in  the  reference  [7]  .  Figure  4  shows  the  conceptual  design  of  the  system  in  which  the  modules 
are  placed  on  the  worksheet  and  then  linked  in  the  appropriate  configuration. 

In  the  scenario  depicted  in  fig  4,  a  natural  guide  star  is  used  to  sense  and  correct  the  tip-tilt  modes 
of  the  atmosphere  whilst  a  laser  beacon  is  used  to  correct  the  higher  order  aberrations  using  a 
Shack-Hartmann  sensor.  In  this  particular  study,  two  equal  turbulent  layers  at  0  and  10  km 
altitude  with  4m/s  and  8m/s  wind  speed  perpendicular  to  each  other,  were  defined  within  the 


8 


ATM  module.  The  wavefront  reconstruction  which  is  performed  after  sensing  the  LGS  with  the 
Shack-Hartmann  sensor  was  carried  out  taking  into  account  only  the  Zernike  modes  ranging 
between  4  and  24  (within  the  REC  module).  A  pure  integrator  was  used  in  the  TFL  module 
(F(s)=  G/s),  the  total  gains  of  the  open-loop  transfer  function  at  0.1  Hz  being  95  Hz  and  159  Hz 
in  the  high  order  and  tip-tilt  adaptive  loops  respectively. 

Notice  that  the  GPR  module  appears  three  times:  once  for  the  NGS  propagation  down  to  the  8m 
observing  telescope,  once  for  the  laser  beam  upward  propagation  from  the  50cm  projector  to  the 
Na-Iayer,  and  once  for  the  LGS  downward  propagation  to  the  8  m  observing  telescope.  Doing  so. 
the  8  m  observing  telescope  parameters  have  to  be  defined  twice,  as  well  as  the  parameters  for 
the  wavefront  correctors  (TTM  and  DMI).  Notice  the  special  module  after  the  TFL  module 
corresponding  to  the  point  where  the  loops  are  closed. 

The  correction  was  considered  in  the  K-band  (X  =  2200  nm  .  AA  =  400  nm )  and  the  PSF 
computed  both  for  the  uncorrected  and  corrected  cases.  Figure  5  shows  the  resulting  long- 
exposure  PSFs,  together  with  the  diffraction-limited  PSF.  One  can  deduce  from  these  plots  that 
the  obtained  Strehl  ratio  is  approximately  0.36.  In  this  particular  case,  it  is  possible  to  form  a 
theoretical  estimate  of  the  compensated  Strehl  ratio  based  on  the  contribution  to  the  overall  phase 
vaiiance  (see  [8]).  The  value  obtained  in  this  way  is  S=0.34  agreeing  closely  with  value  obtained 
from  the  simulation. 


Figure  5.  longitudinal  cuts  through  the  normalized  PSF  (uncorrected,  corrected,  and  diffraction 
limited). 
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2.  BASIS  OF  SIMULATION  MODULES 

In  this  section  we  give  a  brief  description  of  the  model  underpinning  the  major  modules 
contained  in  CAOS. 

Wavefront  Generation 

Three  modules  are  sufficient  to  simulate  the  entire  process  of  wavefront  generation  as  received  at 
the  pupil  of  the  observing  telescope.  Implicit  in  this  procedure  is  specification  of  the 
characteristics  of  the  source  (SRC  module),  the  turbulent  atmospheric  model  (ATM  module),  and 
the  propagation  from  the  source  to  the  observing  telescope  through  the  atmosphere  (GPR 
module). 


2.1.  Atmosphere  Building  (ATM)  Module 

The  ATM  (atmosphere  building)  module  generates  the  turbulent  atmospheric  model  and 
produces  the  distorted  wavefronts  through  the  module  GPR  (geometrical  propagation).  ATM 

therefore  has  no  input  beyond  the  parameters  selected  when  the  parameters  setting  GUI  is 
invoked. 


The  output  of  the  ATM  is  composed  of  a  number  of  turbulent  layers  together  with  their  altitudes. 
Each  turbulent  layer  can  be  physically  simulated  as  a  random  phase  screen  whose  phase  power 
spectrum  has  the  von  Karman  or  the  Kolmogorov  behavior.  If  a  temporal  evolution  is  selected, 
its  base-time  is  fixed  within  ATM  and  then  rules  all  subsequent  simulation  branches.  In  this  case, 
each  turbulent  layer  is  also  affected  by  a  velocity  vector. 

The  first  step  in  building  the  turbulent  atmosphere  is  to  generate  the  phase  screens  that  will 
simulate  the  behavior  of  each  turbulent  layer.  Two  methods  are  currently  available  for  selection: 
a  method  based  on  fast-Fourier-transformation-  (FFT-)  with  boosting  of  low-spatial-frequency 
boosting  [9]  and  an  approach  based  on  random  addition  of  Zernike-polynomials  [10].  Details  are 
given  in  the  associated  references. 

Once  the  phase  screens  have  been  generated,  the  turbulent  atmosphere  is  built  by  arranging  the 
required  number  of  phase  screens/turbulent  layers  with  their  altitudes  and  taking  into  account  the 
Cn  ( h )  profile  chosen  by  the  user.  If  temporal  evolution  is  needed,  this  is  done  once  only.  ATM 

will  then  just  shift  each  of  the  layers  taking  into  account  the  base-time  (minimum 
atmosphere/turbulence  evolution  time)  and  their  associated  velocity  vectors.  If  temporal 
evolution  is  not  selected  ,  no  base-time  and  no  velocity  vectors  are  requested  and  consequently 
each  time  that  ATM  is  called,  it  outputs  a  statistically  independent  ensemble  of  turbulent  layers. 


2.2.  Source  Building  (SRC)  Module 

A  natuial  object  or  an  artificial  laser  beacon  (LGS)  can  be  generated  and  in  both  cases  they  can 
be  defined  as  either  point-like  or  as  an  extended  bidimensional  object.  For  a  natural  object,  the 
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angular  coordinates  of  the  source  (off-axis  and  position  angle),  source  V-magnitude  (background 
magnitudes  for  any  subsequent  noise  computations)  and  spectral  type  can  be  chosen  via  the 
corresponding  parameter  setting  GUI.  The  flux  from  the  natural  source  and  the  sky  background 
are  then  derived  for  the  whole  range  of  optical  Johnson  bands  (from  U  to  M).  The  angular 
coordinates  are  passed  to  the  GPR  module  with  respect  to  a  so-called  main  telescope  having 
coordinates  [0,0,0].  The  actual  telescope  pointing  to  a  given  source  may  be  distinct  from  this.  For 
a  LGS,  its  finite  distance  to  the  telescope  must  also  be  provided  and  the  source  flux  (computed 
from  its  equivalent  V-magnitude)  is  non-zero  only  in  a  narrow  Na-band 
(A  =  589  nm  ,  AA  =  10  nm) . 

2.3.  Geometrical  Propagation  (GPR)  Module 

GPR  peifoims  geometrical  propagation  through  the  atmosphere  by  linear  superposition  of  the 
fluctuations  produced  by  each  turbulent  layer.  The  module  can  be  used  to  simulate  either 
downwaid  piopagation  through  the  atmosphere  (together  with  SRC  and  ATM)  or  propagation 
upwards  through  the  atmosphere.  In  both  cases,  the  simulation  takes  into  account  the  physical 
coordinates  of  the  two  objects  from  which  and  to  which  the  propagation  has  to  be  performed  in 
order  to  establish  which  part  of  each  atmospheric  layer  has  to  be  used. 

For  natural  objects  located  at  infinity,  all  the  wavefront  portions  chosen  at  different  altitudes  have 
the  same  size  as  the  telescope  pupil.  Thus  the  light  beam  propagation  shape  is  cylindrical.  On  the 
other  hand,  a  sodium  LGS  is  situated  at  a  finite  distance  from  the  observing  telescope  and  the 
light  beam  propagation  shape  is  conical.  The  relevant  surface  for  the  upper  layers  is  thus  smaller 
than  for  the  lower  ones.  This  cone  effect  is  simulated  using  a  rescaling  method.  In  other  words, 
before  carrying  out  the  linear  superposition,  each  portion  of  the  wavefront  is  magnified  by  a 

factor  —  where  H Nu  is  the  altitude  of  the  sodium  layer,  and  /;.  is  the  altitude  of  a  given 

turbulent  layer.  The  number  of  photons  from  the  natural  source  and  the  sky  background  are 
computed  using  the  value  of  the  diameter  and  the  obscuration  ratio  of  the  observing  telescope,  as 
well  as  the  respective  fluxes  from  SRC.  The  only  parameters  needed  before  running  this  module 
are  the  coordinates  of  the  observing  telescope  with  respect  to  the  point  [0.0.0],  its  diameter  and 
the  obscuration  ratio. 

2.4  Tip-Tilt  Sensing  and  Correction  (TTS,  TCE  and  TTM)  Modules 

The  TTS  (tip-tilt  sensor)  simulates  the  process  of  imaging  the  wavefront  coming  out  of  the  GPR 
module.  The  designer  can  define  the  number  of  pixels  in  the  detector,  the  field  of  view  subtended 
by  each  pixel,  the  algorithm  to  estimate  the  tip-tilt,  the  sensitive  bandwidth  of  the  detector  and 
the  integration  and  delay  times.  The  module  can  also  take  into  account  different  noise  sources, 
namely  Poisson  photon  and  dark-current  noises  and  Gaussian  read-out  noise.  The  intensity 
pattern  in  the  focal  plane  is  computed  by  using  standard  Fourier-based  scalar  diffraction  theory. 

In  the  case  of  tip-tilt  sensing  with  a  LGS  (generated  with  the  LAS  and  NLS  modules)  each  sub- 
layei  is  tieated  as  an  extended  bidimensional  source,  so  that  for  each  sub-layer  we  compute  a 
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PSF  abei  rated  by  atmospheric  turbulence  plus  the  defocus  aberration  caused  by  the  fact  that  the 
LGS  is  only  focussed  in  a  given  plane. 

In  a  conventional  AO  set-up,  the  TTS  module  is  followed  by  a  tip-tilt  centroiding  calculus 
(TCE)  module  that  takes  the  intensity  pattern  on  the  detector  as  input  and  estimates  the  tip-tilt 
according  to  the  detector  and  algorithm  chosen.  CAOS  currently  provides  barycentre  calculus  for 
a  CCD,  or  quadrant  differential  measurement  for  a  Quad-cell  detector. 

The  required  inputs  into  the  tip-tilt  mirror  (TTM)  module  are  a  GPR-propagated  wavefront 
and  a  set  of  measurements  from  the  TCE  module.  Sometimes,  the  designer  may  wish  to  first 
filter  these  in  the  time  domain  using  the  TFL  module.  A  tilted  plane  is  then  subtracted  from  the 
wavefront  in  order  to  correct  the  measured  tip-tilt.  No  dynamic  behaviour  of  the  tip-tilt  mirror  is 
explicitly  considered  although  this  is  possible  with  a  suitable  choice  of  time  filter  (see  next 
section). 

2.5  Time  Filtering  (TFL)  Module 

The  time  filtering  (TFL)  module  implements  a  recursive  digital  filter  in  the  time  domain. 

Because  most  AO  systems  work  in  closed  loop,  the  TFL  module  can  be  used  to  emulate  a  servo- 
control  law  and  time  filter  the  instantaneous  estimate  of  the  residual  wavefront  requiring 
correction  before  applying  it  to  the  wavefront  corrector.  The  design  of  the  servo-control  law  is 
critical  in  order  to  ensure  the  closed-loop  stability  and  provide  optimal  performance  for  a  given 
AO  system,  atmospheric  conditions  and  reference  source.20'21  The  TFL  module,  however  Js  not 
limited  to  simulation  of  a  control  law  and  can  be  used  whenever  time  filtering  is  needed.  For 
example,  noise  filtering  in  open  loop  configuration  can  be  simulated  in  this  way. 

The  equivalent  Laplace  domain  transfer  function  for  the  digital  filter  design  can  be  selected  from 
any  of  the  following:  (i)  a  generic  analog  filter  expressed  in  terms  of  gain,  zeros  and  poles,  (ii)  a 
single  pole  at  zero  frequency  with  a  user  defined  gain  G,  i.e.  a  pure  integrator  with  H(s)  =  G/s. 
and  (iii)  a  proportional-integrator-derivative  (PID)  filter  given  by- 

«(«)-*, +£-+*, -!i i- 

S  S  +  COq 

wheie  K,,,K,  and  Ktl  are  the  gains  for  the  proportional,  integral  and  derivative 

components  of  the  filter,  the  low-pass  correction  term  -(°n  with  cut-off  frequency 

tu0+.v 

is  introduced  to  filter  out  the  contribution  of  the  high-frequency  noise  in  the  derivative 
component.  In  each  case,  the  Z-transform  of  the  corresponding  digital  filter  is  computed  via  use 
of  the  bilinear  approximation,  and  the  recurrence  relationship  is  shown  within  the  parameter 
setting  GUI.  One  can  also  plot  the  Bode  diagram  of  the  equivalent  analog  filter.  The  TFL  module 
code  accepts  a  command  type  structure  data  as  input,  in  which  the  present  values  of  the  discrete¬ 
time  signals  to  filter  are  stored.  The  same  digital  filter  is  applied  to  each  signal,  and  the  result  is 
passed  in  another  command  type  structure  data  as  output. 
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2.6  Shack-Hartmann  Sensing  (SHS)  Module 

The  SHS  module  computes  the  image  of  the  guide  object  (i.e.  natural,  point-like  guide  star, 
extended  laser  beacon  etc)  as  viewed  by  each  subaperture  of  the  Shack-Hartmann  lenslet  array.  It 
then  combines  them  into  the  full  sensor  image  resized  to  the  CCD  scale,  integrates  the  images 
over  the  integration  time  and  adds  the  noise.  Its  input  is  the  wavefront  and  its  output  is  the  CCD 
sensor  image  integrated  over  the  required  time. 

The  Hartmann  lenslet  geometries  may  be  selected  from  the  following - 

i )  A  square  geometry. 

ii)  Hexagonal  geometry. 

iii)  Radial  geometry  with  equal  area  subapertures. 

i  v)  A  generalised  radial  geometry  in  which  the  number  of  rings  in  the  array  and  the 

number  of  subapertures  per  ring  may  be  specified 
v)  Customised  geometry  defined  by  means  of  a  user-supplied  file  (which  must  contain 
the  sensor  geometry  array  and  optical  axis  positions  of  the  lenslet  in  this  array). 

The  default  setting  is  to  assume  that  the  pupil  is  centered  on  the  sensor.  However,  the  size  of  the 
pupil  lelative  to  the  sensor  and  its  relative  shift  can  be  specified  if  required. 

CCD  camera  characteristics  can  also  be  specified  within  this  module.  The  number  of  camera 
pixels  corresponding  to  one  subaperture,  the  angular  size  of  a  CCD  pixel,  the  subaperture  angular 
size  (if  there  is  a  field  stop  placed  in  front  of  them),  the  integration  time  and  read-out  time  and 
the  working  wavelength  band  can  be  set.  It  is  also  possible  to  consider  the  noises  (in  the  same 
way  as  in  the  TTS  module),  and  to  use  a  standard  Gaussian  calibration  fiber  or  a  specific  one 
(adjustable  FWHM  and  shape).  Differential  elongation  of  LGS  tridimensional  spot  and  Rayleigh 
scatteiing  cone  on  each  subaperture  have  also  been  implemented. 

An  example  of  a  simulated  CCD  sensor  image  is  shown  in  Fig.  6  for  a  square  7x7  Shack- 
Hartmann  sensor,  with  8x8  pixels  under  each  subaperture.  The  aberrated  wavefront  corresponds 
to  third  order  astigmatism. 

2.7  Centroiding  calculus  module  (CEN) 

The  CEN  module  takes  the  output  from  a  SHS  CCD  image  and  defines  the  pixels  over  which  the 
estimate  of  the  subaperture  centroid  is  to  calculated.  The  calculated  centroid  positions  can  be 
linearly  related  to  the  local  slopes  of  the  wavefront  averaged  over  the  subapertures.  One  can 
choose  to  rebin  the  CCD  camera  pixels  over  which  the  centroiding  is  made  and  can  set  the  rebin 
f  actor  to  any  submultiple  of  the  linear  number  of  pixels  per  subaperture.  One  can  also  decide  to 
use  the  fibre  reference  measurements  as  zeros  for  CEN  output  measurements  or  to  give  the  raw 
(non-i eferenced)  data.  The  extracted  pixels,  on  which  centroiding  is  computed,  can  be  centred  on 
the  subapeiture  optical  axis  (usual  case)  or  on  the  reference  measurement  central  positions 
(useful  in  case  of  large  static  aberrations). 
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Figure  6.  Shack-Hartmann  wavefront  sensor  (left),  representation  of  the  spots  (centre)  and  of  the 
measured  slopes  (right)  for  a  square  7x7  subaperture  sensor. 


2.8  Reconstruction  (REC)  Module 

Input  to  the  REC  module  is  the  centroiding  measurements  which  are  output  by  CEN.  The  REC 
module  can  be  employed  to  perform  three  basic  functions  - 

(i)  To  visualize  the  measured  approximation  of  the  wavefront  for  the  Fried  geometry.  In 
this  case,  the  output  is  the  local  wavefront  phase  reconstructed  at  the  actuator 
positions. 

(ii)  To  compute  the  mirror  commands  to  be  sent  to  the  time  filtering  module.  In  this  case, 
the  output  is  the  mirror  command  voltages. 

(iii)  To  evaluate  the  modal  coefficients  of  the  measured  wavefront.  In  this  case,  the  output 
of  the  module  is  the  modal  coefficients  of  the  specified  basis  functions.  Any  of  the 
three  bases:  mirror  modes,  Zernike  circular  polynomials  or  a  Karhuenen-Loeve  basis 
may  be  specified.  The  interaction  matrix  can  be  automatically  calibrated  by  the 
simulation  or  can  be  provided  by  the  user.  Command-to-mode  and  mode-to-command 
transfer  matrices  are  also  computed  from  the  mirror  influence  functions. 


2.9  Deformable  Mirror  (DMI)  Module 

The  DMI  module  builds  the  deformable  mirror  geometry  as  well  as  its  influence  functions,  and 
uses  the  commands  computed  by  the  reconstructor  and  time  filtering  modules  in  closed  loop  to 
correct  the  incident  wavefront.  It  therefore  requires  two  inputs:  the  incident  wavefront  to  be 
corrected  and  the  commands  coming  from  the  reconstruction  (REC)  or  time  filtering  (TFL) 
module.  DMI  has  two  outputs:  the  corrected  wavefront  and  the  mirror  geometry  parameters 
needed  by  the  reconstructor  module  to  build  the  interaction  and  passage  matrices  (number  of 
active  actuators  and  influence  functions  of  the  mirror).  One  can  specify  the  geometry  of  the 
mirror,  its  influence  function  characteristics,  and  its  dynamic  behavior.  If  required,  one  can  also 
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specify  the  size  and  central  shift  of  the  pupil  on  the  mirror,  the  angle  of  incidence  of  the 
wavefront  and  the  mirror  azimuth  angle  in  its  own  plane  relative  to  the  plane  of  incidence. 

Currently,  2  pre-defined  piezo-stack  mirror  geometries  are  provided.  These  are  square  and 
rectangular  geometries  with  the  number  of  actuators  and  spacing  along  x-  and  y-axes  specified 
by  the  user.  A  custom  geometry  can  also  be  defined.  The  dynamical  behavior  of  the  actuators  is 
characterized  by  the  gain  (displacement  induced  by  a  unit  voltage  command),  by  the  maximum 
admissible  actuator  stroke  and  by  the  time  delay  between  reception  of  the  command  and 
application  of  the  stroke.  Hysteresis  is  not  currently  modelled. 

2.10  Data  analysis-system  performance 

Whilst  the  aim  of  CAOS  is  not  to  provide  extensive  data  analysis  tools,  basic  system 
performance  is  provided  through  the  following  modules  - 

PSF  -  a  point  spread  function/image  computation  module 
DIS  -  a  generic  data  display  module 
SAV  -  for  saving  data  generated  during  the  simulation 
STF  -  a  structure  function  computation  module 

The  specific  functionality  of  these  modules  can  be  found  in  the  user  documentation. 

2.11  CAOS  -  CONCLUSION 

Our  experience  to  date  has  shown  that  CAOS  is  a  powerful  tool  for  the  simulation  of  adaptive 
optics  systems  enabling  studies  of  scientific  issues  to  be  carried  out  with  much  less  effort  than 
would  be  required  by  writing  one’s  own  code. 

Particularly  attractive  features  are  i)  There  is  no  need  for  the  user  to  get  involved  in  low  level 
coding  li)  The  well-defined  and  clear  use  of  modules  enabling  simple  modeLling  and  design  iii) 
In  those  cases  where  CAOS  cannot  provide  the  required  functionality,  the  user  may  contribute 
their  own  module  provided  it  follows  the  established  protocol  (though  this  has  not  been  carried 
out  in  our  case).  Disadvantages  of  the  CAOS  approach  in  our  experience  are  i)  The  overall 
flexibility  means  that  the  generated  code  is  not  optimized.  Execution  time  can  be  long  and  it  is 
demanding  on  computer  memory  resources  (in  general,  128  Mbytes  Ram  should  be  considered  a 
minimum),  ii)  The  software  is  not  as  stable  as  would  be  ideal. 

In  conclusion,  however,  we  believe  that  CAOS  ranks  as  an  excellent  tool  for  AO  studies.  It’s 
continued  development  and  improvement  and  more  widespread  use  should  mean  that  it  will 
become  an  increasingly  important  resource  for  the  AO  community. 
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3.  Modal  wavefront  projection  functions  with  minimum  error  norm 

In  this  following  section,  we  present  an  article  submitted  for  publication  on  June  8lh  2000  to 
JOSA  A.  The  paper  is  a  theoretical  contribution  which  shows  that  we  may  derive  vector 
polynomials  to  effectively  restore  orthogonality  to  the  representation  of  wavefronts  as  a  modal 
expansion  when  the  wavefront  sensing  technique  provides  measurement  of  the  local  gradients  of 
the  wavefront.  A  variational  formulation  shows  that  the  polynomials  are  optimum  as  they  will 
guarantee  that  the  error  propagation  in  the  wavefront  reconstruction  is  minimum.  We  anticipate 
that  this  subject  may  be  of  relevance  to  many  workers  in  the  fields  of  wavefront  sensing,  adaptive 
optics  and  optical  testing. 


VARIATIONAL  SOLUTION  FOR  MODAL  WAVEFRONT 
PROJECTION  FUNCTIONS  OF  MINIMUM  ERROR  NORM 
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ABSTRACT 

Common  wavefront  sensors  such  as  the  Hartmann  or  curvature  sensor  provide 

measurements  of  the  local  gradient  or  Laplacian  of  the  wavefront.  Expression  of  wavefronts 
in  terms  of  a  set  of  orthogonal  basis  functions  thus  generally  leads  to  a  linear  wavefront 
estimation  problem  in  which  modal  cross-coupling' occurs.  Auxiliary  vector  functions  may 
be  derived  which  effectively  restore  the  orthogonality  of  the  problem  and  enable  the  modes 
of  a  wavefront  to  be  independently  and  directly  projected  from  slope  measurements.  Using 
variational  methods,  we  derive  the  necessary  and  sufficient  condition  for  these  auxiliary 
vector  functions  to  have  minimum  error  norm.  For  the  specific  case  of  a  slope-based  sensor 
and  a  basis  set  comprising  the  Zernike  circular  polynomials,  these  functions  are  precisely  the 


Gavrielides  functions. 


18 


1.  INTRODUCTION 

In  problems  of  wavefront  estimation,  it  is  common  to  expand  the  wavefront  in  terms  of  a  set 
of  orthogonal  functions  or  modes  - 

*00=2V.(*)  (in 

which  obey  the  orthogonality  relation  - 

J/j,(x)P<(x)^^  =  5m  (i.2) 

D 

and  D  denotes  the  domain  of  integration. 

The  most  commonly  used  set  of  functions  in  optics  is  the  Zernike  circular  polynomials. 
Speaking  generally,  orthogonal  bases  are  desirable  because  they  allow  the  modal 
coefficients  to  be  evaluated  by  simple  integration  of  a  product  of  two  functions  over  the 
domain,  D  using  the  orthogonality  relation  given  by  eq.  1 .2.  However,  since  Hartmann 
sensors  (and  shearing  interferometers)  provide  estimates  of  the  gradient  of  the  phase  rather 
than  the  phase  itself,  the  appropriate  model  is  — 

0.3) 

where  the  coordinates  \  .  denote  the  position  of  the  slope  measurements  within  the  pupil.  In 

this  case,  the  orthogonality  of  the  wavefront  basis  cannot  be  exploited.  Modal  cross¬ 
coupling  will  occur  and  the  coefficients  must  be  obtained  by  solving  an  inverse 
problem/overdetermined  system  of  linear  equations.  Many  workers  have  examined  this 
problem  in  a  search  for  optimal  solutions  [1-6]. 

The  evaluation  of  modes  by  direct  integration  could  be  restored  if  we  can  derive  a  set  of 
auxiliary  vector  functions  F,  (x)  which  are  orthogonal  to  the  gradients  of  the  basis  functions. 


To  establish  the  necessary  conditions  for  such  a  set  to  exist,  consider  that  the  divergence  of 
these  vector  functions  gives  the  basis  functions  i.e.  we  have  the  relation  - 
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V-F,(x)  =  />(x)  (1.4) 

Using  the  orthogonality  relation  eq.  1 .2  in  eq.  1.1  and  substituting  for  Pj  (x)  in  eq.  1 .4.  the 
modal  coefficients  are  given  by  - 

a;=J<p(x)V-F((x)r/2.v  (1.5) 

D 

This  may  be  written  as  - 

fl/=J  V’{F-  (x)<p(x)}r/2.v-  f  V(o(x)-F/  (x)  d2.x  (1.6) 

o  D 

Applying  the  divergence  theorem,  we  obtain  - 

=  J  <p(x)  F(.  (x)-  c/I  -  J  V<p(x)-F;  (x)  d2x  (1.7) 

D  D 

Clearly  if  the  required  set  of  functions  F,  (x)  satisfies  the  relations  - 
V.F;(x)  =  /^(x) 

F,(x)-«/l  =  0  (1.8) 

where  F;  (x)-  d  1  denotes  the  normal  component  of  the  function  to  the  closed  contour  C  of 
the  integration  domain  D,  we  then  have  - 

«/  =  -  J  V<p(x)-F;(x)  d2x  (1.9) 

D 


In  other  words,  we  may  evaluate  the  modes  by  direct  integration  as  required.  It  has  been 
shown  that  there  are,  in  fact,  a  number  of  possible  sets  of  vector  polynomials  which  may  be 
used  in  eq.  1 .9  [7-9]  and  there  is  no  unique  solution  to  the  problem  described  by  eq.  1 .8. 


2.  VARIATIONAL  FORMULATION 


In  the  remainder  of  this  paper,  we  will  address  the  question  of  how  to  derive  the  optimum  set 
of  vector  functions  Fi  (x)  .  By  optimum,  we  mean  that  set  of  vector  functions  which  will  give 
the  minimum  error  or  noise  propagation  in  our  estimate  of  the  wavefront  given  by  eq.  1.9. 

The  estimate  of  the  kth  modal  coefficient  of  the  wavefront  will  be  given  in  practice  by  - 

at =  Jd(x)F*(x)  d2x  (1.10) 

D 

where  d(x)  =  V<p(x)  +  v(x)  and  v(x)  is  the  additive  noise  vector  at  x. 

If  we  may  assume  a  noise  process  having  zero  mean  and  covariance  cr2S  (x  -x')  (which  is 

reasonable  for  Hartmann-type  sensors),  it  is  straightforward  to  show  that  the  variance  in 
radians"  associated  with  our  estimate  of  the  k  mode  is 

V(F>^(x)|V,  (1.11) 

D 

The  noise  propagator  is  thus  defined  as  N  (Fk )  =  J|F(  (x)|2  d\x  and  depends  only  on  the 

D 

particular  choice  of  the  vector  functions  through  the  volume  defined  by  |Fa  (x)|2  in  D.  The 
ensemble-averaged  mean-square  error  associated  with  the  estimator  is  then  - 


2 

where  the  latter  term  simply  corresponds  to  unestimated  modes.  As  the  conditions 
specified  in  eq.  1.8  do  not  determine  a  unique  solution  for  F*  (x),  there  are  many 
possible  solutions  for  each  mode,  each  having  a  different  associated  noise  propagator. 

The  key  point  is  that  we  seek  to  minimise  the  noise  propagators  expressed  by  eq.  1.10 
subject  to  the  constraint  equation  and  boundary  condition  expressed  by  eq.  1 .8.  This  is 
a  variational  problem.  Accordingly,  we  introduce  the  Lagrange  multiplier  function 
K  (x)and  seek  to  minimise  an  objective  function  given  by  - 

Q  =  J|F*  (x)f  d2x+ j K  (x){V  •  F*  (x)-/>  (x)} d2x  (1.12) 

D  D 

We  now  take  first  variations  in  Q  with  respect  to  F(  (x)and  set  these  to  zero  to  find 
the  stationary  points  of  Q  - 

e  =  2j5F1(X).F,(x)</!.r+jA,(x)V.5F,(x)</!A-  =  0  (1.13) 

D  D 

Since  V  {A,  (x)Ft  (x)}= A,  (x)V  F<  (x)+Ft  (x)- VA(  (x),  this  may  be  expressed  as  - 
2  =  2j«F,  (x)-F,  (x)AtJv.[  A.  (x)5F,  (x)}c/=.r- Ja,  (x)aF,  (x)efl=0 

D  D  D 

(1.14) 

Applying  the' divergence  theorem  to  the  second  integral  in  eq.  1.14  and  rearranging 
we  obtain  - 

G  =  J{ 2F,  (x)- VA,  (x)}.<SF,  (xy’ar+JSF,  (x)  A,  (x)Jh  -  0 

D  c 

(1.15) 

The  boundary  condition  Ft  (x)-dl  =0  from  eq.  1 .8  ensures  that  the  last  integral  in  eq. 

1 . 15  is  identically  zero  and  we  arrive  at  the  simple  condition  for  the  FA.  (x)  to  have 
minimum  error  norm  - 
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F*  (x)  =  VA,  (x)  (1.16) 

In  other  words,  the  F<  (x)  are  an  irrotational  set  of  functions  -  i.e.  must  be  expressible  as 

the  gradient  of  an  associated  scalar  function.  Note  we  have  dropped  the  factor  of  2  as  this 
will  have  no  effect  on  the  functional  form.  Substituting  eq.  1.16  into  the  constraint  given 
by  eq  1.8  we  find  that  the  associated  scalar  field  is  determined  by  the  Poisson  equation 
with  Neumann  boundary  condition  - 

V2A,(x)  =  Pa  (x)  xeD 

VA,(r)-n(r)  =  0  reC  (U7) 

This  is  the  main  result  of  our  analysis.  The  Poisson  equation  with  Neumann  boundary’ 
conditions  has  a  unique  solution  and  the  minimum  error  norm,  modal  projector 
functions  can  be  obtained  by  solving  eq.  1.17  for  the  given  set  of  basis  functions, 

Pk  (x).  Such  a  set  of  functions  was,  in  fact,  explicitly  derived  in  analytic  form  for  the 

Zernike  circular  polynomial  basis  by  Gavrielides  in  1982  [10]  although  the  noise 
propagation  behaviour  of  these  functions  was  not  addressed  in  his  original  work.  For 
completeness,  their  analytic  form  is  given  here  - 

The  Gavrielides' vector  polynomials  are  defined  in  polar  coordinates  as  follows  - 
F,  M)  =  < (r,0)  f +  Fj  ( r,G )  t?  where 
For  m  *  0 


F;  =^2(n  +  \)Tj(r)cos(m&) 

if  i  is  even 

Fr  =~yl^(n  +  \)  Tj"  (r)  sin(/m9) 

if  i  is  odd 

F$  =~V2(n  +  0  Qn(r)  [-/«sin(mt?)] 

if  i  is  even 

Fg  =—^2(n  + 1)  Qjj  (r)  [mcos(;m?)] 

if  i  is  odd 

23 


n-m 

i 


C  (0=71 

4  ,v=0 


C(s)(/i-2.v  +  2) 


/z  +  m 


•s  + 1 


V  n  —  m 


—  j+l 


n-m 

I 


e:w^S 

tf/zd  where 

cw= 


A  - 

cm 


n  +  m 


—  s  + 1 


V 


n  -  in 


~  J  +  1 


[r'"_l  -r"-2l+l] 


[(« -2i-  + 2) -mr"-2t+l] 


A 


-I'  (/!  -5)! 


5! 


'  n  —  m  ^ 

f  n  +  m  ^ 

s 

2  J 

s 

2  J 

For  zrz  =  0,  vv>£  have 


F;(r)=l^T\  £ 

7T  ,v=0 

r;(0=o 


-l'(/t-.v)!  r"-2v+l 


Yn  > 

si 

—  s 

f 

? 

- 

A  ~  J 

_ 

(U8) 


3.  CONCLUSION 

Given  the  inherent  simplicity  of  modal  projection  from  slope  measurements  accordina 
to  eq.  1.7,  one  may  speculate  that  these  functions  should  see  more  widespread  use  in 
wavefront  sensing  and  estimation  problems.  We  anticipate  that  the  knowledge  of  their 
minimum  noise  propagation  behaviour  will  make  them  a  simple  and  attractive 

approach  to  wavefront  estimation  from  slope  measurements  and  lead  to  their 
increasing  use. 

It  is  interesting  to  note  that  the  problem  we  have  discussed  of  obtaining  a  set  of 
auxiliary  vector  functions  which  gives  minimum  wavefront  error  norm  has  a  close 
analogy  with  minimum  physical  principles.  Eq.  1 . 1 7  may  in  fact  be  identically  used  to 
describe  the  equilibrium  shape  of  a  membrane  (with  fixed  boundary  described  by  the 
closed  contour  C)  over  which  a  force  Pk  (x)acts  to  produce  vertical  displacement 
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K  (x)-  In  this  case,  J|Fa.  (x)| 2  d2x  represents  the  potential  energy  stored  in  the 

D 

membrane  and  will  be  equal  to  the  work  done  in  bringing  the  membrane  to  it's 
equilibrium  position  J P  (x) Xk  (x)d2x  only  in  the  case  of  conservative  fields  -  i.e. 

D 

when  Fk  (x)  =  VAJt(x)  [II]. 
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